tfd(ind) = c_d0(ind)+g_d0(ind)+exports0(ind);
tfs(ind) = l0(ind)+sum(ind2,iof0(ind2,ind));

parameter
aa(ind,ind2)
;

aa(ind,ind2) = iod0(ind,ind2);

scalar mm; mm = 0;
scalar colflg; colflg = 0;
scalar rowflg; rowflg = 0;
scalar mflg; mflg = 0;
scalar iflg; iflg = 0;


parameter
chk_c(ind)
sum_c(ind)
chk_r(ind)
sum_r(ind)
scal_c(ind)
scal_r(ind)
;

while( (mm<500) and (mflg < 1),
mm = mm+1;

* COLUMNS
iflg = 0;

chk_c(ind2) = tfd(ind2) + sum(ind,aa(ind2,ind));
sum_c(ind2) = tfs(ind2) + sum(ind,aa(ind,ind2));

scal_c(ind2) = sqrt(chk_c(ind2)/sum_c(ind2));

loop(ind2,
if (  (abs(sum_c(ind2)-chk_c(ind2)) >= 0.00000000001),
iflg = iflg + 1;
);
);

loop(ind,
aa(ind,ind2) = scal_c(ind2)*aa(ind,ind2);
);

if( (iflg = 0),
colflg = 1;
);

iflg = 0;
chk_r(ind) = tfs(ind) + sum(ind2,aa(ind2,ind));
sum_r(ind) = tfd(ind) + sum(ind2,aa(ind,ind2));
scal_r(ind) = sqrt( chk_r(ind)/sum_r(ind) );

loop(ind,
if (  (abs(sum_r(ind)-chk_r(ind)) >= 0.0000001),
iflg = iflg + 1;
);
);

loop(ind2,
aa(ind,ind2) = scal_r(ind)*aa(ind,ind2);
);

if( (iflg = 0),
rowflg = 1;
);

* Check for convergence

if( (rowflg = 1) and (colflg = 1),
mflg = 1;
);


);


display
iflg
scal_c
scal_r
colflg
rowflg
mm
mflg
;

ex(ind) = (tfd(ind) +  sum(ind2,aa(ind,ind2)))  - (tfs(ind)+sum(ind2,aa(ind2,ind)));

iod0(ind,ind2) = aa(ind,ind2);
